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ABSTRACT 

An approximation involving cubic spline functions for 
parameter estimation problems in the Euler-Bernoulli-beam 
equation (phrased as an optimization problem with respect to the 
parameters) is described and convergence is proved. The 
resulting algorithm was implemented and several of the test 
examples are documented. It is observed that the use of penalty 
terms in the cost functional can improve the rate of convergence. 


Research supported TrJ part by the "Fonds zur Forderung der 
wissenschaf tlichen Forschung", Austria, No. 4534 and in part by 
the National Aeronautics and Space Administration under NASA 
Contract No. NASI— 17070 while the second author was a visitor at 
the Institute for Computer Applications in Science and 
Engineering, NASA Langley Research Center, Hampton, VA 23665. 


“T ' 


i 




1. Introduction: Statement of the Problem and Notation. 


In this investigation we are concerned with the parameter 
dependent equation 


(1.1) 

y tt = - ^(q^y) - D 2 (q 2 D 2 y t ) 

- q 3 y t + q^y , 

(1.2) 

y(0, •) = y Q , y t (0,*) = y 1 , 


(1.3) 

boundary condition , 


where 

y = y(t,x) , 0 < x < 1 and where D 

stands for 


differentiation with respect to x. Together with appropriate 
boundary conditions, equations (1.1) - (1.2) arise in modelling 
the tranverse vibrations of a thin and elastic beam, where struc- 
tural (q^ * 0) or viscous (q^ * 0) dampling may be present [6, 
Chapter 17) . In an effort to provide a mathematical description 
of the motion of large flexible systems, continuum models of the 
type (1.1) - (1.3) are receiving anew much attention in the 
present literature [2, 3," 7, 10, 18]. 

Different boundary conditions generally require a somewhat 
different analysis and we therefore restrict ourselves to the 
following boundary conditions of type k, k=l or 2 or 3 : 


(k=l): y(t,0) = y xx (t,0) = y(t,l) = y xx (t,l) = 0 

(simply supported beam), 

(k=2): y(t,0) = y x (t,0) = y(t,l) = y x (t,l) = 0 

(beam clamped at both ends). 
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(k= 3 ): y(t,0) = y x (t,0) 


for t > 0 . 


y xx (t > 0) “ ^xxx^’^ " ® 

(cantilevered beam). 


Let us assume that (1.1), ( 1 . 2 ) together with boundary conditions 
of type k are the mathematical model of a physical system, for 
which observations $(t..,x.) , i=l,...,l and j=l,...,m , 

are available and which are expected to correspond to the 
solution y(t i ,x j .;q) of (1.1) - (1.3) evaluated at (t i ,x j .) , 
with the parameter 


q = ^ q l ,q 2 ,q 3 ,q 4 ’ y o ,y l^ 

chosen "correctly". Since it is mathematically untractable and 
practicably not feasable to require £(t.,x.) = y(t.,x.;q) for 

and we formulate the following optimization 

problem 

(P) minimize £ |$(t.,x.) - y(t- ,x.jq )| 2 over 

i,j J J 

q e Q subject to y satisfying (1.1) - (1.3) , 

and refer to (P) as parameter estimation problem. Here y(*,*;q) 
denotes the solution of (1.1) - (1.3) in a sense that will be 
specified below and Q, the set of admissable parameters is 
given by 

Q C Q 1 x Q 2 x l" * L°° x H* j x H° , 
where Q 1 = ^ ^ H 2 : q^x) j a } , with a > 0 , 





Q 2 = {q 2 6 H 2 : q 2 (x) 5 0} 

H (l)^ q l) = {v e H 2 : v(0) = v(l) = 0} , 

2 2 

H( 2 )(q^) = {v e H : v(0) = v(l) = v'(0) = v'(l) = 0} , 
H 2 5) (q 1 ) = (v G H 2 : v(0) = v'(0) = 0} . 


The set Q is given a topology by taking 


,ci: “weak - "weak - L* * L* * H 2 * H° . 


2 p 

Here H we ak denotes the space H endowed with the weak topology. 


The parameter estimation problem (P) is an infinite dimensional 
problem, and its approximation by a sequence of finite dimensional 
problems that are numerically implementable is the focus of this 
paper. Here we study approximation of (P) by projection on finite 

dimensional subspaces. By factoring the operators D 2 (q 1 D 2 ) and 

2 2. * 

D (q 2 D ) into second order operators, the smoothness requirement 

on the functions in the subspaces yrill be that of H 2 -regularity . 

In tlqis respect our treatment differs from earlier work (see 

[2,3,7]), where it is assumed that the functions in the subspaces 

are sufficiently smooth, so that they lie in the domain of the 

generator of the semigroup associated with (1,1) - (1.3), in 

particular that they are H^-functions . In special cases this 

requirement can be circumvented in [2,3*7] by transforming the 

equation into a higher dimensional system of equations. The present 

approach has the advantage that we can give a complete convergence 


r 
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theory for the approximation of (P), which includes, for example, 
cubic spline approximations. 

We also study, theoretically and numerically, the effect of 
adding a penalty term to the quadratic fit-to-data criterion 
involved in defining (P). The use of such a penalty term seems 
to be justified for practical problems and proved to be numerically 
useful, since it provides some robustness in the search for a 
minimum of (P). In some examples the computing time was reduced 
and in others less care was required to tune the parameters in 
the minimization routine when a penalty term was added. This will 
be illustrated with an example at the end of the paper. 

The "observations" for our test problems were obtained by a 
generalized Crank-Nicolson-scheme , exept in very special cases, 
when the solution could be calculated analytically. Since numerical 
approximations for the trajectory of (1.1) - (1.3) do not seem 
to be readily available in the literature, we also compared the 
Crank-Nicolson-scheme and a cubic spline-projection algorithm 
for approximation of the trajectory with an example where the 
analytical solution is available. It turns out that the projection 
algorithm is superior to the Crank-Nicolson-scheme, both with 
respect to time and accuracy. 

In our model equation we have not included a forcing function 
f(t,x). The results of this paper can easily be extended to in- 
clude such a term, with techniques that are fairly routine and 
are discussed in detail in [5] for second order equations. 
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In section 2 we formulate (1.1) - (1.3) as an abstract differential 
equation in a Hilbert space and prove a general convergence 
result for projection schemes that is uniform in the parameters. 
This parameter dependent state convergence result is used in 
section 3 to study approximation of (P). Section 4 is devoted 
to a short discussion of the use of penalization techniques 
that proved to be numerically very useful for the present problem. 
We have carried out extensive tests for a cubic spline scheme 
and an example is discussed in section 5 • A detailed description 
of the scheme and of the phenomena that were observed will appear 
elsewhere . 

The notation that we employ is rather standard and we only make 
a few comments. By |* | and (•,•) we denote the usual norm and 
inner product in L . Other norms are denoted by a subscript and 
ll *!lx use ^ to specify the operator norm in the space X. 

We will not distinguish between column and row vectors. As usual, 
a(A) and p(A) stand for spectrum and resolvent set of a 
closed operator A and A| D denotes the restriction of A to the 
set D. 


T 
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2. Abstract Formulation and Parameter Dependent State Convergence. 


In the first part of this section we restrict our attention to 
the equation 

(2.1) y tfc = - D 2 (q 1 D 2 y) - D 2 (q 2 D 2 y fc ) , 

y(o, • ) = y Q , y t (o, • ) = y 1 , 

boundary conditions of type k , 


where 0 £ x £ 1, t ;> 0, y Q e y^ G H° and k=l or 2 or 3- 

The general case given in (1.1) - (1.3) will follow from the 
results concerning (2.1) by adding a suitably defined bounded 
perturbation. At first we also assume y Q and y^ to 
be fixed. We reformulate (2.1) in a function space setting and 
give the convergence result that is the essential tool for the 

approximation of the parameter estimation problem (P). Several 

. . 2 
preliminaries are required. For q^ e the spaces ^^(q^) 

that were introduced in the previous section are endowed with the 

inner product 


<V y W> 


(q 1 D 2 v j D 2 w) 


2 2 
In this way become Hilbert spaces and the norm lv| q = 

2 2 

<v,v>q is equivalent to the usual H -norm on H^) > 

moreover this equivalence is uniform as q^ varies in bounded 
subsets of Q ; for details see [ 7 , Theorem 1.8), for example. 
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Further we introduce the sets 

hJ 1} = {v e H 4 : V(0) = v"(0) = v(l) = v’’(l) = 0} , 

hJ 2) = {v 6 H 4 : v(0) = v'(0) = v( 1 ) = v'(l) = 0} , 

hJ 3) = (v G H 4 : v(0) = v'(0) = v"(l) = v"'(l) = 0} . 

Let Y be the space H° endowed with the weighted inner product 

-1 

(v,w) y = (q 1 v,w) , 


with G and define closed linear operators C(q^) from 

H° to Y by 

p 

dom C( ) = ^(k) 5 
C(q 1 )v = q a D 2 v . 


Note that 


<v j w>q = (C(q 1 )v, C(q 1 )w) Y for v,w e H (k) 


and that 


( 2 . 2 ) 


I C ( q. ) v l v s k min q.(x)lv! 0 , 

1 1 xG[0,l] 1 H 2 


for all v G ^(k) an d a constant k independent of q^ and v. 


Lemma 2.1. The adjoint 
H° is given by 





of 


C(q a ) 


mapping from Y to 
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dom C* 1) (q 1 ) = {v e H 2 : v(0) = v(l) = 0} , 
dom C* 2) (q 1 ) = H 2 , 

dom C*^^) = {v'6 H 2 : v(l) = v'(l) = 0} , 

and 

^(k) ^ q l = ® v * 

Proof . This technical result follows, with minor modifications 
due to the weighted H°-norm of Y, from the calculations in [11, 
pg. 169-171]. 


Lemma 2.2. Let = C (k) ^ q l ) • Then A (k)^ q l^ is 

selfadjoint and strictly positive in H° and dom A ^ k ^(q 1 ) = 


is a core for C(q^) . 


Proof . By von Neumann's theorem [11, pg. 275] 5 A (k)^ q l^ 
selfadjoint and dom A^jCq^) is a core for CCq^) . The 
remaining assertions can easily be checked. 


Corollary 2.1, 


and of H 


00 


The inclusions of H 
in H° are dense. 


00 


in H 


00 


(in the H 


00 ' 


Proof . Recall that by the definition of a core, the set 
{(VjCCq^Jv) : v e H^ k ^> is dense in the graph of CCq^ . For 
each (v,C(q^)v) with v £ H^ k ) there exists a sequence 

v n e H (k) w ^ th v n * v and C ^ q l^ v n * in Si nce 

2 i| 

the H^ k ^-norm and the graph norm of C(q^) are equivalent, H^ k ^ 


norm) 


y 


2 2 
is dense in H^ k j . Finally, H (k) 

the claim follows. 


is clearly dense in h and 


Next we define for qj € Qj a family of operators A(q 2 ) by 

Z) 

dom A(q 2 ) = and 


A( k )( q 2 ) v = D 2 (q 2 D 2 v) 


To write (2.1) as a system of first order equations, we introduce 
the Hilbert spaces 


X = H 2 k) ( qi ) x H° , 


with the inner product denoted by ( • , • . In the notation of X 

as well as with other symbols we now adopt the convention that 

contents permitting the indices k and q may be dropped . For q G Q 

1| /] 

let be gi ven by domA^^Cq) = H^ k j x H^ k j and define 


■*(k) (q) 


- ( ° 


(k)^ q l^ 


-A 


1 ) • 

(k) (q 2 } / 


We shall need the fol-lowing technical lemma. 


Lemma 2.3. For each q^ G there exists a constant K such 

that 

I A(q 1 ) I a Klvl a , 
x H 

. If, moreover, q^ varies in a bounded (with 


q 

for all v G H^ k ) 
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2 

respect to the H -norm) subset of then K is independent of q^. 

Remark 2.1, In this .lemma the requirement q^ e h 2 is needed 
in an essential way. 


Proof of Lemma 2.3. The proof of such an inequality is fairly 

standard and we shall therefore only sketch it. Since 

!A(q 1 )v! £ K 1 I v I for some constant K. depending only on Iq | 

**• 1 

it suffices to show that lACq^vl S K 2 ID 4 vl , where K 2 satisfies 
the specified properties. We have 


a 2 |D^v| 2 £ iq^vi 2 = |A(q 1 )v| 2 - |D 2 q 1 D 2 v + 2Dq 1 D 3 v,| 2 
1 

2 | (q 1 D^v ) (D 2 q^D 2 v + 2Dq 1 D 3 v)dx . 


For any e > 0 we find 


1 = 2 


iq 1 D iJ vi I D 2 q^D 2 v + 2Dq 1 D 3 vldx 


' 1 

S 2lq l' H l^ J ll)i,v| lD 2 q 1 D 2 v|dx + 2 l Dq^ | lD 3 v|dx) 


lq ! 1 2 

^ 2lq | . ( 2e I D v I 2 + |D 2 v| 2 + 

H 1 4 E L°° 


lq l |2 l 

— ID 3 v I 2 ) 
2e 


Using a basic estimate [17, pg.l8] we have for any 0 < e 1 < l 
and an appropriately defined constant K independend of q^ : 


A. . 2 


I S 2 I q 1 1 ^ [ 2e I D H v | ^ + Elq 1 l 2 e~ 1 (z^ 2 \v\ 2 + e^D^vl 2 )] 
^ H 


A.. 2. 


IT 
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In a similar way we find that 

II = lD 2 q 1 D 2 v + 2Dq 1 D^vl 2 S I q 1 1 2 2 &( e^ 2 I v I + e^D^vl) , 

with £ independent of q^. Those two estimates are used in (2.3) 
and for appropriately chosen constants e and we conclude 
that lA(q^)vl £ KjlD v| . This ends the proof. 

Proposition 2.1. For q € Q the operator ^(^(Q) is densely 
defined and dissipative. Its closure ^(q) is maximal dissipative 
and generates a C Q -contraction semigroup T(t;q) on X. 


Proof . By Corollary 2.1 ^ ^(q) is densely defined. 

Dissipativity of , A'(j c )( ( l) is easily verified. We note that 

ii h 

for real X > 0 and x sufficiently small, (X-«A^ k > (q x H^^) 

h ii ii li 

is dense in H (j c ) x ^(k) * In f° r a11 ( a > b ) 6 ^(k) X H (k) we 

find that (X-A^j (q) ) (^) = (^) is equivalent to v = Xu-a and 


(2.4) 


(X^+A(q 1 ))u = -XA(q 2 )u + f 


with f = Xa + A(q 2 )a + b G H° . Let (X 2 +A(q^))u = w . Then (2.4) 
becomes 

(2.5) w = -\A(q 2 )(X 2 +A( qi )) _1 w + f , 


or w = -xA(q 2 )A(q^) ^A(q^ ) (x 2 +A(q^ ) ) + f . We solve for w e H° 

in* (2.5). Note the A(q^ ) (x 2 +A(q^ ) ) -1 is a bounded operator in H° 
uniformly in \ > 0, and by Lemma 2.3 A(q 2 )A(q^) is a bounded 

operator as well. Consequently a fixed point argument guarantees 


T 
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the existence of a solution w e H° of (2.5) for some \ 


2 —1 ti 

sufficiently small. We put u = (X^Afq^)) w e H (k) 


and the 


required density of the range of (\ Q -J^ k j(q)) is shown. 


Now the dissipative operator iA^ k ^(q) admits a closure ^( k )(q) 
[12, pg. 86]. Moreover the range of (x o -^ k ) (q) ) is necessarily 
closed in X [12, pg. 86]; but it is also dense and so it coincides 
with X. This implies that the closure jv^ k ^(q) of A-( k )(q) is 
a maximal dissipative operator and that it generates the semi- 
group T(t;q) on X. 


Remark 2.2. If q ? = 0 , then dom j^ k j(q) = and 

the action of ■A'( k )(q) coincides with that of A- k (q) • Moreover 
A"( k )(q) is skew adjoint in this case, as can be seen by employing 
von Neumann’s theorem [11], for example. - For arbitrary q 6 Q 
it is simple to verify that dom A-^Cq) c H 2 k ) x H 2 k ) and that 


< A' (k )(q)(v ) ^b )> X = <V ’ a> qi " (Q!D 2 U,D 2 b) - (q 2 D 2 v,D 2 b) 


2 2 

holds for every (a,b) e H ( k ) x H (k) * 


Recall that for z 6 dom A^ k ^ we have 


( 2 . 6 ) 


TCtjqJz = ^( k ) (q)T(t;q)z , for t > 0 


This equation is formally equivalent to (2.1) with T( t ;q) (y Q ,y ) = 

(y,y t ) • 
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We proceed to describe an operator theoretic form of a Galerkin 
approximation to (2.1). Let X N be a sequence of finite 
dimensional subspaces of and denote by X^ (X^) the 

space X N endowed with the H°-(H^ k) (q 1 )-) topology. Note that X N 
depends on k as well. Further we introduce the orthogonal 
projections pJJ : H° * xjj and : H^ k) (q 1 ) - xjj . We put 
X N = X^ x X^ and define : X - X N by 


0 

o : P N 


Whenever the dependence of P^ and P N on q is relevant 
write P^(Q 1 ) and P N (q) . 


For q 1 6 let 


C N (q 1 ) = C(q ± ) I X N 


Clearly C^(q^) G i(X^,Y) and the adjoint C^(q^)* satisfies 


(C N (q 1 )u,v) y = (u,C N (q^ )*v) Q , for all u G X^ , v G Y . 

H 


We next define 


A N ( qi ) = .C N (q 1 )*C N (q 1 ) 


and note that A N (q 1 ) G L(X N ) and 


(2.7) (A N (q^ )u,v ) = (C N (q 1 )u,C N (q 1 )v) Y = (C(q 1 )u > C(q 1 )v) y , 


for all UjV G X . Here we defined the approximating operators 
on X N only, but there is an obvious extension to all 
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of H° given by Ag Xt (q a ) = (C(q 1 )P^) # C(q 1 )P^ . For q 2 € Q 2 
define operators A^(q 2 ) : X^-X^ by 

(A N (q 2 )u,v) = (q 2 D 2 u,D 2 v) for all u,v e X N . 


Finally we introduce 
*oc> (q) ; - 

the operator 

^( k) (q) e £(x n ) 

approximating 

■*(k) (q) * 

/ ° 

^ ~ A Oc) ^ q l ^ 

/ )• 
A (k) (q 2" 



The Galerkin approximation to (2.1) can now be introduced as 

<A'^j t )(q)z (t;q) , 

P N (y 0 ,y 1 ) . 

Note that (2.8) is an equation in the finite dimensional space X N . 

Lemma 2.k. The operators *0^ k )(q) , q 6 Q , generate contraction 
semigroups T N (t;q) on X N for each N. We have T N ( t ;q)P N (y Q ,y ) = 
z N (t;q) . 

Proof. From the construction of vA^ k )(q) it follows that 
(^ k )(q)z,z) x s 0 for all z e X N . This and the finite dimen- 
sionality of X N imply the result, see [12, pg. 90]. 


( 2 . 8 ) 


<nr = 


z N (0;q) = 
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Theorem 2.1. Let and q° G Q , N=l,2..., and assume that 
pJJ * I strongly in H° and p!J(q°) - I strongly in H^(q°) . 
If moreover (qj,qlj) * (qj.q^) weakly in H 2 x H 2 , then 

T N (t;q N )(P N (q N )z » T(t jq°)z in X , 

uniformly in t as t varies in bounded subsets of [0,<») and for 
every z G X . 

The approximation of the general equation (1.1) - (1.3) will 
follow trivially from this theorem, see Corollary 2.1 below. 
Before we give the proof of Theorem 2.1 we establish some 
technical lemmas. 


Remark 2.3. Strictly speaking the above convergence result is 
not stated precisely, since the spaces in which the operators are 
defined vary topologically. A more precise statement would be 
given by 

IT N (t;q N )J> N (q N )3 N z - tr N T(t;q°)zl M - 0 as N - - 

X(q N ) 

uniformly in t in bounded intervals of [0,»). Here Cf N *. X(q°) - X(q N ) 

is the canonical isomorphism. Since the set {iq^i : N=l,2,...} 
is bounded, the spaces X(q^) are topologically equivalent 
and no confusion should arise by the abbreviated statement of 
the theorem. 


Lemma 2.5. For 
that of A( k) (q) 


q € Q , the resolvent set of 
contains 0. 


iA-(k)(q) 


as well as 
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Proof. Let (u,v) € X N such that A N (q)(u,v) = (0,0) . Then 


<* N (q)(£),(y)> x = 0 for all (x,y) e X N 


This is equivalent to <v,x> - <u,y> - (A N (q ? )v,y) = 0 . First 

q l q l 

let y = 0 , and conclude that v = 0 . Next putting x = 0 with y 
arbitrary implies u = 0 and in particular 0 6 p(Ji N (q)) . 


Similarily, since 0 e pCACq^)) it is simple to show that 0 is not 
an eigenvalue of JKq) and that the range of JKq) includes 

i| 2j 

H(k) x H(k) and is therefore dense in X. We will show that 

1 Zj 4 

*A(q) I H(k) x H (k) con tinuous (in X). Let (x,y) e 

H (k) x H (k) • Then O = is equivalent to 


(2.9) v = x and u = -A(q 1 ) -1 (y+A(q 2 )x) . 


We have 


( 2 . lo ) lA(q 1 )“ 1 ylg = (C(q 1 )A(q 1 ) _1 y J C(q 1 )A(q 1 ) 1 y) Y = (A(q a ) V^y) S 


const (q ± ) ly I . 


Using [11, pg. 169] we find in a similar manner 


(2.11) lA(q 1 ) 1 A(q 2 )xl = lC(q 1 )A(q 1 )" 1 A(q 2 )x l y = lC(q 1 )C(q 1 )~ 1 (C*(q 1 )) 1 A(q 2 )xl Y 


= l(C*(q 1 )) _1 D 2 (q 2 D 2 x)l y = 1-^- q 2 D 2 xl S const (q)ixl q ^ . 
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The estimates (2.10), (2.11) together with (2.9) imply 
l(u,v)ijj s Kl(x,y)i x , for an appropriately defined constant 
depending on q and indepent of (x,y) and the claim is verified. 


Lemma 2.6. For q 6 Q and z e we have j (q) _1 z 
P N A' (k) (q)" 1 z . 

Proof. Let z = (x,y) and define ) = iA N (q) -1 z and 

M 

( ) = <A(q) z . We need to show that 


( 2 . 12 ) 


u N = P?u and v N = P N v 
1 o 


For every N we have 

N 


.N 


-f-A N ( q l)u N- A N( q2)v N j =(y) , 


and in particular v N = x . Next we use Remark 2.2 and find that 

2 2 

for every (a,b) £ x 


(2.13) 


<S(q)(“),(§)> = <v > a> qi " (Q^ u,D 2 b) - (q 2 D v,D 2 b) = <x,a> q ^ + (y,b) . 


This implies v = x , so that the second equality in (2.12) is 
satisfied. Since -0 £ p(«A(q)) the first equality in (2.12) will 
be verified, once we show that P^u is the (unique) solution of 
A N (q^)u N = -A N (q 2 )x - y in X N . In the following calculation 
we use (2.13) with a = 0 and b € X 1 ^: 


T 
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-(A N (q )P^u,b) = -<P?u,b> = -<u 3 b> = (y,b) + (q 2 D 2 x 3 D 2 b) = (y,b) + (A N (q 2 )x 3 b) 

J. i. -L 4^ 4^ 

Since b £ X N is arbitrary this equality implies the result. 


The convergence statement in the following lemma has to be inter- 
preted in the sense explained in Remark 2.3* 


Lemma 2.7. Let q N , q° € Q with q N - q° weakly in Q. Then 

'^( k )^)" 1 - ^( k )^ 0 )” 1 in x • 

h h 

Proof . Let (x 3 y) £ ^(ic) x ^(k) anc * define 

(*> = and ($) = Mq N ) _1 (y) . 

Then v = v N = x , u N = -A(q^) - 1 (A(q 2 )x+y ) and u = -A(q°) (A(q 2 )x+y) . 

From Lemma 2.3 we conclude that {lu^l j. : N=l,2,...} is bounded. 

H 

Therefore there exists a subsequence of u^ again denoted by u N , 

4 2 4 

converging weakly in H and strongly in H to an element z € ^(k) * 

2 

We have for every w £ ^(k) 

-(q^D 2 u N ,D 2 w) - (q^D 2 x,D 2 w) = (y,w) . 

Taking the limit in this last inequality we get 

-(q°D 2 z 3 D 2 w) - (q 2 D 2 x,D 2 w) = (y,w) . 

This implies z = u . The usual subsequence argument can be used 

N 2 

to show that the original sequence u w ■* u in H ( ^ ^ ( q^ ) . Since 
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{ll ( A'" 1 (q N )ll 


M : N=l,2,...} is uniformly bounded in N the 

X(q N ) 


h 2i 

result follows from density of H^ k ^ x H^ k ^ in X. 

For the ease of the reader we state a version of the Trotter-Kato 
approximation theorem that is readily applicable for the proof 
of Theorem 2.1. 


Lemma 2.8. [19] Let q N G Q , q° G Q and assume that 

lP N (q N )zl .. - | z I for every z G X . Then 

X(q N ) X(q°) 


(2.14) lp N (q N )T(t;q°)z - T N (t ; q N )J> N (q N ) z I - 0 

X(q N ) 


for- every z 6 X if and only if there exists a X 6 n p(A-^^(q^)) n 
p(«A( k ) (q° ) ) such that 


(2.15) l(X-^ k) (q N )) _1 P N (q N )z - F N (q N )(x-^ r (q 0 )) - ^! -0 

(k) X(q N ) 


for every z G X and 


(2.16) 


MT N (t;q N )ll xj £ Me ut 
X(q N ) 


for positive constants M and o independent of N. 

Again, a remark analogous to Remark 2.3 has to be made; (2.14) 
for example, should read 


- T N (t;q N )? N (q N )3 N zl -0 

X(q N ) 


lP N (q N )3 N T(t;q°)z 
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Proof of Theorem 2.1. First we note that weak convergence of q^ 
N o 

implies that q™ is bounded in H and therefore the I • | „-norm 
1 _N 

2 

is uniformly in N equivalent to the H -norm. Since 



|P ^(q^)v-v! N £ lP^(q°)v-v 

' N 

S K IP^(q°)v-vl 


q l 

q l 

q l 

for some constant independent 

of 

p 

N and v 6 , it follows 

that 

lP l( q l ) v-v 1 N * 0 as N * 
q l 

oo 

and 

(2.17) 

lP N (q N )z-zl - 0 

as 

N ~ for every z £ X 


X(q ) 


Moreover, for 


v 6 H 


2 

00 


we find 



lvl nI = l ,P K )vl 


N 



and this last term converges to 0 as N - «* . 

Therefore l? N (q N )zl „ - Izl ^ for every z G X . To verify 

X0q N ) X(q°) 

the convergence result we use the estimate 


lT N (t;q N )? N (q N )z - T(t;q°)zl x S 

(2.18) 

IT N (t;q N )f> N (q N )z - ? N T( t ;q° ) z I x + l? N T(t;q°)z - T(t;q°)zl x 

The set (T(t;q°) : t € [0 ,t]} is compact for every x > 0 and 
therefore (2.17) implies that the second term on the right hand 
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side of ( 2 . 1 8 ) converges to zero uniformly in t on bounded 
intervals. To show convergence of the first term on the right 
hand side of (2.18) we use Lemma 2.8. In view of the previous 
calculations it suffices to verify (2.15) with \ = 0 . By 
Lemma 2 . 5 and 2.6 we have for z 6 X 


.A N (q N )-V*(q N )z „ * 

X(q N ) 

(2.19) lT N (q N )»A'(q N )‘ 1 r N (q N )z - T N (q N )A(q N ) _1 zl + 

X(q N ) 

lP N (q N U(q N ) _1 z - P N (q N )>(q 0 ) _1 zl 

X(q N ) 


This last term converges to 0 by Lemma 2.7. Since 
{lip (q^)Mq^) ^11 : N=l,2,...} is bounded, (2.17) implies conver- 
gence to 0 of the first term on the right hand side of ( 2 . 19 ) 
and the result is proved. 

We now return to equation (1.1) - ( 1 . 3 ) and define operator an 
£( k) (q) by dom A| k) (q) =dom£ (k) (q) and 


£( k )U) = & (k) (q) + B(q ) 


where 


B(q) 


.-q- 


0 


“B 

Since is a bounded perturbation of ^ (k )(q) it generates 


a semigroup, denoted by T B (t;q) for every q e Q . Further we 
define approximating operators ^£®(q) on X^ by 



22 


(2.20) *(£®(q) = *(k) (q) + > 

so that the Galerkin approximations to (1.1) - (1.3) become 


( 2 . 21 ) 


^z N (t;q) = ^^®(q)z N (t;q) , 
z N (0;q) = F N (y ojyi ) . 


The solution to (2.21) is given by the semigroup Tg(t;q) 
generated by A^£)(q) . 


Corollary 2.1. Let q N , 

and assume that - I 

o 

in H (k)(^i) • Then 

T^(t;q N )? N (q N )(y 


q° £ Q with q N - 
strongly in H° and 


N y N ) - 
o ,y l J 


T B ( t ; q° ) (y 


q° weakly in Q, 
P^(q°) » I strongly 


,y 1 ) in X , 


uniformly in t in bounded subsets of [O,oo). 


Proof , To verify the claim we employ the triangle inequality: 


lT^(t;q N )? N (q N )(yN }y N) _ T B ( t ;q° ) (y Q ,y ± ) 


' X S 


lT£(t;q N )? N (q N )(yJJ,y!J) - T g( t ;q N )? N (q N ) (y Q .y^l ^ + 
lTN(t;q N )P N (q N )(y o , yi ) - ^(t ;q° ) (y Q ^ ) l x . 


The operator norms IIJ , ^(q^)B(q^)|| 


v, are uniformly bounded in N. 

X(q N ) 
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Consequently (see e.g. [ 15> Theorem 3-1*1]) llTg(t ;q^)P^(q^) II 
is bounded uniformly in N and t in bounded intervals of [0,~). 
This implies convergence to zero of the first term on the right 
hand side of (2.22). As for the second, we note that for z e X 

lT N (q N )B(q N )z - B N (q°)P N (q N )z| „ £ lB(q N )z - B(q°)f }N (q N )z I N 

X(q N ) X(q N ) 

S lB(q N )z - B(q°)z I + lB(q°)z - B(q°)? N (q N )z I I 

X(q N ) X(q N ) 

and these two terms converge to zero by assumption. A simple 
perturbation result as stated in [19] implies the convergence 
of the last term in (2.22) and the corollary is verified. 

Remark 2.*J. In [ 10] the author applies the Trotter-Kato theorem 
to demonstrate the convergence of Galerkin and finite element 
approximations to second order hyperbolic systems. Our analysis 
here is similar, but more complicated due to the parameter 
dependence of the operators (and spaces) and the perturbation 

p p 

term D (qgD y^.), which is of maximal order. The main step in 
the proof of Theorem 2.1, the convergence of the resolvents, is 
reduced to convergence of the projections and the results of 
Lemma 2.6 and 2.7 - 

Remark 2 . 5 » Let us briefly compare the convergence result of this 
paper to the results of [2,3>7]. The main difference is given by 
the fact that Banks and Crowley in their work take 
the finite dimensional subspaces as subsets of the domain of the 
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generator. This generally requires more smoothness of the 
elements in the subspace than our approach does and in addition 
more boundary conditions have to be met, (compare h! to H? ) . 

\ *W V K ) 

However the original operator ^ k) (q) can be factored so that 

this smoothness requirement can be relaxed, see [2,7]. In that 

case parameter dependent convergence was only shown for 

q 2 = q 3 = 0 in 17]. Furthermore in their application of the 

Trotter-Kato theorem. Banks and Crowley show convergence of the 

approximating generators on dense subsets of X, which implies 

convergence of the resolvents, whereas we show convergence of 

the resolvents directly. This allows for the generality that 
N N 

q a , q 2 need only to converge weakly as opposed to strong con- 
vergence of , q 2 required in [3]. As a consequence, the 

usual compactness assumption on Q in parameter estimation problems 
(see (HI) below) is fairly weak in the present case. 

Remar k 2.6 . We have chosen three specific boundary conditions 
here. Other "natural" boundary conditions, as for example 
combining a simply supported beam with a clamped end, can be 
treated by the technique that we used here, provided that it is 
guaranteed that HCvfl defines a norm on dom C which is equivalent 
to the H —norm and that the result of Lemma 2.3 holds. For a 
useful discussion of various boundary conditions see [8] , for 
example. 
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3. Approximation of the Parameter Estimation Problem. 


We return to the optimization problem (P) of Section 1. It will 
be convenient to introduce the notation 

(u(t;q) j v(t;q) ) = T B (t;q)z , (u N (t;q), v N (t;q)) = Tg(t;q)J> N z , 

for z = (y y ) e X . The abstract formulation of (P) is given by 

p 

(AP) minimize Y i£(t.,x.) - u(t.,x.;q)l as q varies in Q , 

i,j J J 

and a sequence of optimization problems which will be shown to 
approximate (AP) is defined by 

(AP)^ minimize Y i#(t.,x.) - u N ( t . } x . ;q ) i B as q varies in Q . 

i,j 1 J 1 J 

As in Section 2 we take X N ^ H^ k) and X N = x!J x xJJ . 

The following hypotheses will be used: 


(HI) Q is a compact subset of Q (in the norm of 
H weak * H weak - L~ * I*" * H^ k) x H°) , 


(H2) Pq - I strongly in H° and P^q^) * 1 strongly in 
*2 


H (k) for every q.^ contained in the first component of Q. 


For the fit-to-data criteria in (AP) and (AP) N we also use the 
notation 


T 
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and 


J(q) = - u(t.,x.;q)i 2 

jN (q) = ~ u N (t i ,x. ;q) | 2 


Lemma 3 • 1 . Let (HI) hold. Then (4P; has a solution q^ € Q for 
each N. 


Proof. Prom (2.21) and Gronwall's lemma it follows that q n « q° 
in Q implies T 3 (t,q n )z - T^’^(t;q°)z in X and consequently 
u ( t ;q n ) - u N (t;q°) in C ( uniformly in bounded intervals of t. 

The result then follows from the special form of the fit-to-data 
criterion and (HI ) . 

Below N p and M p denote subsequences of N that tend to » as p - » 

Theorem 3.1. Assume that (HI) and (H2) hold and let q* e Q be a 

• • . -N _N 

limit point of (q } , so that lim q p = q* is Q . Then q* is a 

solution of (AP) and 


(3.1) 

and . 

(3.2) 


N 


N 


>n,- n p - ] 


N 


T ^(t^q P )P (yjf* ,y^) - T(t;q*)(y*,y*) uniformly in 

bounded intervals of [0,»), 


N N. 

J P (q“ P ) - J( q* ) . 


Pr oof . Compactness of Q implies the existence of a subsequence 
with the specified properties. For every p we have 


N n N r, N 

J P (q P ) S J P (q) , 


(3.3) 
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for every q 6 Q . We use Corollary 2.1 and take the limit on 
both sides of (3.3). This implies J(q*) £ J(q) for all q 6 Q . 

The remaining assertions follow trivially from Corollary 2.1 . 

Remark 3.1. It is easily seen that a more general fit to data 
criterion could be employed as well. In fact, let 
J(q ) = I(u( • jq) ,v( • ;q) ,£) and J N (q) = I(u N ( • ;q) ,v N ( • ;q) ,$) 
where I( • ,$) :C(0,T;X) - R is continuous and y denotes some 
appropriate measurements. Then the results of Theorem 3-1 remain 
unchanged . 

Remark 3.2. The numerically inclined reader might ask why we 
did not attempt to demonstrate any rate of convergence, provided 
the projections converge to the identity with a certain rate. 

Certainly the methods of section 2 allow to demonstrate a rate 
of convergence of the approximating semigroups to the semigroup 
T (t;q) on X, uniformly in q 6 Q . It seems to be difficult, 
however, to obtain any estimate on the rate of convergence of 
the solutions q^ of (AP)^ to q*. Such a rate was only proved 
in a special case of an elliptic equation so far 
[9], The technique used in [9] does not generalize to other types 
of equations and moreover it requires that every solution of the 
infinite dimensional optimization problem satisfies a certain 
regularity property.* 

Remark 3.3. If the set Q in Theorem 3.1 is finite dimensional 

* T 

K. Kuni'sch and L. W. White have studied this regularity assumption 
recently in a forthcoming paper "Regularity of the diffusion coefficient 
m parameter estimation problems.'- 


T 



(e.g. a set of cubic polynomials with unknown coefficients as 
parameters), then (AP)^ is a finite dimensional problem. If Q 
itself is infinite dimensional, then a further approximation 
of the elements in Q is required. We shall not pursue this 
problem since the techniaues and difficulties that arise 
do not differ essentially from those for parabolic problems of 
second order which have been investigated in some detail in [4,13]. 
For a result on the approximation of variable- coefficients in 
beam equations see [ 3 ]. 


4. The penalized fit-to-data criterion. 

The set Q in Theorem 3.1 is constrained in two distinct ways. 
First there are pointwise constraints given in which 

guarantee wellposedness of the equation (1.1) r- (1.3) and 
secondly there are constraints inherent in the compactness 
assumption that give existence of solutions of (AP) and (AP) N . 
In our numerical implementation of the problems (AP) N , however, 
we use unconstrained optimization algorithms because they are 
more simple to handle, - and because they work. To parallel at 
least partially what is done in practice we briefly discuss 
the use of a simple penalization method. Another aspect calls' 
for the use of a penalization method as well: In practical 
examples it is likely that a priori knowledge of the system 
that is modelled leads us to expect the unknown parameter to 


29 


lie within certain bounds and this knowledge should be used in 
the fit-to-data criterion. In all the numerical tests that we 
carried out a deliberate use of' a penalty term decreased the 
cost of the calculations. For a discussion of the penalty method 

for parameter estimation problems in second order parabolic 
equations we refer to [14]. 

To simplify the presentation we assume the parameters (q 
y y ) to be known and restrict our attention to estimating q^. 

We agree to drop the index "l” henceforth so that q = q 1 • For y > a 
we introduce the set 

(4.1) Q = (q G H 2 : Iq I ? S y } . 

1 H 

Note that n is nonempty closed and convex and therefore 

a weakly closed subset of H^. Since n is also bounded, 

this set is weakly compact in n and Theorem 3-1 is applicable, 
provided (H2) holds. We need the following hypothesis: 

(H3) The functional ij> : H - jR is weakly lower semicontinuous 

* A 

and satisfies \f;(q) =0 if and only if q G and 

ip(q) 00 if iql p 00 

In addition to (AP) and (AP) N we introduce the problems 

(APP) M minimize J(q) + Mtjj(q) as q varies in Q 1 , 

and 

(APP) M,N minimize J N (q) + Mip(q) as q varies in Q 1 . 


T 
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Lemma 4 . 1 . Let (H3) hold. Then (APP)^ and (APP)^’^ have solutions 

q M and qjj in Q r The set { iqjjl 2 : N=l, . . . ; M=l, . . . } is bounded. 

H 

Proof • Since is radially unbounded, the existence of q^ and q„ 
trivially follows. Moreover for every q € n we have 

0.2) 0 £ J N (qJ{) + M ^(qJ}) S J N (q) . 

But (J N (q) : N=1 , . . . } is bounded by Theorem 3.1 and therefore 

( ^ ( q M ) J N-l , . . . ; M=1 , . . . } is bounded. Using (H3) this implies 
the claim. 

Theor em 4.1. Assume that (H2) and (H3) hold and let q* be a weak 

-N N 

limit point of { q M > so that q M P - q* weakly in H 2 with 

P 

N P " “ 5 M P * “ ' Then q * 6 ®i n Q i is a solution of (AP) and 
N _N N 

T P (t;q M P )J> P (y o ,y a ) - T(t;q*)(y o , yi ) in X , 

uniformly in bounded intervals of [0,~), and 
N N_ N 

J P (q M P ) + M <Kq M P ) - J(q*) as p r - . 

P K p 

Proof - The short proof is quite standard but it is included here 

for the sake of completeness. Prom (4.2) it follows that 

_N N 

lim ^(q^ ) = 0 . Since 0 s iMq*) S lim inf iKq M P ) =0 we have 
* p P-“ M p 

A 

q e Q a n Q 1 . The following inequality holds for each q € ^ n Q 1 


Taking the limit p - - in (4.3) and using Theorem 3-1 we 
arrive at 

(l|.l|) J( q* ) £ J(q) for- every q 6 Q a n Q a . 

This implies that q* is a solution of (AP). The final assertion 

* 

in the statement of the theorem follows from (^.3) with q - q . 


5. An Example. 


In this section we present numerical data for the approximation 

of unknown -coefficients , where the exact solution of the 

N 2 

differential equation is known. The subspaces X c are 

chosen as cubic spline functions, see [16, Chapter 4 ] with the 
basis reduced by two degrees of freedom for k=l or 3 (four 
degrees for k=2) to account for the boundary conditions. By 
using elementary estimates on spline interpolation, it is simple 
to show that (H2) holds for boundary conditions of type k=l,2 or 3 
Once a basis for X N is chosen a matrix representation for the 
ordinary differential equation corresponding to (2.21) is readily 
derived. Depending on the parameters, this ordinary differential 
equation may be stiff or also highly oscillatory. In the first 
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case the NAG-impler.er.tancn of the Gear algorithm was used to 
solve the ordinary differential equation, in the latter we 
either used the Gear algorithm or the NAG-implementation of 
the Adams-Bashford algorithm. The minimization problem (AP)^ 
was solved by the IMSL-implementation of the Levenberg-Marquardt 
algorithm. Some experience has to be obtained in tuning the 
parameters of the Levenverg-Marquardt algorithm as applied to 
the present problem. We had very good success using a three- 
procedure: from one step to the next the convergence 
criteria for the Levenberg-Marquardt algorithm as well as the 
tolerance for the error test determining the variable stepsize 
and order of the Gear (or Adams-Bashford) algorithm are sharpened, 
and the approximated parameters of the previous step are used 
as start up values in the succeeding one. All computations were 
carried out in double precision in FORTRAN on the UNIVAC 1100/81 
computer at the computing center of the Technical University of 
Graz. - Examples similar to the one presented here involving 
quintic and cubic spline approximations are described in [2,3,7]. 
As mentioned in the introduction, cubic spline approximation 
could he used in [ 2 , 3 , 7 ] only after transforming (1.1) - (1.3) 
to a system of three equations. In this case the structure of 
the approximating system of the ordinary differential equations 
is different from ours. 

In our numerical experiments we took synthetic data as our 
measurements S^(t i ,x j . ) . For the example shown below the exact 
value of the analytical solution was used, in all other examples 
the data were generated by a finite difference scheme. 


We now turn to the specific equation 


y tt = -fljA - <32 D ‘ , 5 r t • 

y(0,x) = q.gSinnx , 

:5.d 

y t (0,x) = 0 , 

boundary conditions of type 1 , 

where q^,q 2 an( ^ Qg are taken to be constants. 

The solution of (5.1) for t a 0 and 0 s x s 1 is given by 

(5.2a) u ( t j x ) = exp(-a 1 t)(C^exp(a 2 t ) + C 2 exp(-a 2 t)}sinnx , 

if A = n^q 2 - ^q^ > 0 , 

2 2 

(5.2b) u(t,x) = exp(-a 1 t)(C 1 sin a^t + C 2 cos a 2 t}sinnx , 
if A < 0 , and 

(5.2c) u(t,x) = (C^exp(-a 2 t) + C^t exp(-a 1 t ) )sinnx , 

if A = 0 . Here a^ = (n**q 2 )/2 , a 2 = (n^\HA T )/ 2 , 

cj = qg (a^+a 2 )/ ( 2a 2 ) , C 2 = qg ( a 2 _ a^ )/2a 2 ) , = a^qg/a 2 , 

C 2 = q g = C3 and = qg(a 1 ~l) . 

In the numerical experiments the "correct" value of the parameter 
vector was taken to be 
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Each solution of (AP)^ requires the approximate solution of 

(5.1) via (2.21) many times (see "function evaluations" below). 

It is therefore important that the algorithm for the approximation 

of the partial differential equation is a powerful one. We 

compared the numerical data given by the implementation of (2.21) 

• N 

with X taken as cubic spline functions as mentioned above, 
to a finite difference scheme which is constructed by generalizing 
an explicit scheme discussed in [1, pg. 280]. For (5.1) the 
projection-spline-scheme was superior with regards to accuracy 
as well as computing time. (The experience with all other 
numerical examples showed superiority of the projection - over 
the finite difference scheme as well.) 

In the tables below we show the relative error in percent of 
the observed data, which is defined as 

relative error in percent =1 g- x - act val u e - approximation ! 

/ exact value r uu ' 

Table 1 shows the relative error (in percent) of the solution 
calculated by a finite difference technique in the time interval 
[0,2]. ‘We took an equidistant grid of [0,l] of grid length N -1 , 
and a grid along the t-axis of grid size 2J -1 . The numerical solution is 
given at (t,x) = (2, .5) where the relative error in percent 
E^( 2 .,. 5 ) reached its maximum. 



35 


Table 1. Relative error i.p. of state at (2., .5) (finite difference). 


N/J 

loo 

looo 

2ooo 

l*ooo 

CPU-time for J=l*ooo 

1* 

2.759 

1.716 

1.659 

1.631 

I.07I* sec. 

8 

1.535 

0.1*933 

0. 11366 

0.1*082 

2.563. 

16 

1 .2*49 

0.2071 

0. 15oU 

0. 1221 

5.686 

32 

1.179 

0. 1367 

0.0800 

0.05166 

11.929 

6 k 

1.161 

0.1192 

o.o62l*7 

0.031*11* 

2U .1* 17 

128 

1.157 

0.11U8 

0.05809 

0.02976 

1*9.391* 

Table 2 

gives 

the data 

for the 

same forward problem calculated 


by means of the projection method with cubic spline subspaces, 
using the Gear algorithm to solve the ordinary differential 
equations and with N defined as above. The data are significantly 
better. The rate of convergence seems to be like N _1 * in this 
example . 


Table 2. Relative error i.p. of the state at (2., .5) (cubic splines) . 


N ET (2^.5) maximum of absolute error CPU-time 


1* 

.81*55 * lo 2 

.7885 x lo" 1 * 

1.998 sec. 

8 

.1*959 x lo -3 

.1*616 x 1o~ 5 

u.263 

16 

-1* 

•3o17 x 10 

.2827 x lo" 6 

10.252 

Next we show 

some results on 

the approximation of 

the ’’unknown 


parameter values (q^>§2>§g) • The Levenberg-Marquardt minimization 





algorithm requires a start-up value and it was taken to be 
(0,0,0) for every N, although in practical examples one would 
generally take the solution of (AP)^ ^ as the start-up value 
for (AP) N . 

Table 3 contains the relative error in percent when each of the 
parameters is searched for separately. 


Table 3. 


Relative error i.p. for 


q l’ q 2 * q 3 separately. 

E 16 CPU-time for N=U,l 6 


*2 

q 6 


.17 x 

10 2 

• 

CD 

X 

lo ** 

.27 X 

1o- U 

io.6. 

6 0.7 

.11 X 

1o“ 2 

.5!* X 

lo -1 * 

X 

C\J 

lo" 1 * 

n. 14, 

65.8 

.2h x 

lo- 2 

. 18 X 

10- 3 

.92 X 

lo'" 14 

8.8, 

5 o .1 


For each value of N the total number of iteration for the "three 
step" Levenberg-Marquardt algorithm was 7 for q^, 9 or 10 for q 2 
and 5 for qg. A total of 19 function evaluations was necessary 
for the estimation of q^, 22 for and 15 for qg. When pairs of 
parameter values were searched simultaneously the results were 
not significantly different and the calculations were roughly 
twice as costly. Finally we give the data, when (9} >92 ,q 6^ were 
searched simultaneously. Again the start-up value is (0,0,0). 
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Table 4 . Relative err.or i.p. for (q 1 ,q2>Qg) simultaneously. 

E* E 8 E l6 CPU-time for 

.57 * lo” 1 .31 x lo" 2 .11 x 1o“ 2 N=4: 36.602 

.57 x lo -1 .31 x lo -2 .11 * lo” 2 N=8: 5o.6l*9 

.26 X lo” 14 .1*0 X 1o -U .52 x lo -1 * N=l6 : 137. ool 

Depending on N there were 52-59 function evaluations necessary 
and 16-17 iterations of the Levenberg-Marquardt algorithm. When 
the start-up values in the three parameter search were taken as 
large as t*5,-5»5-) the algorithm was shown to be still 
convergent . 

We also made experiments to study the effect of noise-corrupted 
observations. Table 5 shows the relative error in percent of 
the estimated parameters for N = 4, when Gaussian noise with 
standard deviation o is added to the exact solution ?(t^,Xj) . 
For N=8 or 16 the relative error was not significantly different 
from that for N=4. Moreover , . comparing the relative error of the 
estimated parameters under exact measurements with those of the 
corrupted ones, we observe that for o=.0i the relative error in 
percent increases roughly speaking hundred times and for o=.l 
about thousand times. Rows 1-3 give the relative errors in percent 
when the parameters are estimated separately, rows 4-9, when 
they are searched for in pairs and rows 10-12 show the results 
for simultaneous estimation of (q^,q 2 >qg)> A cross in the table 
indicates that the results were obviously worthless. 


q 6 
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Table 5. Relative error i.p. for N = 4 and noise-corrupted data. 


a = .00 a = . o 1 a = . 1 


*1 

*2 

*6 

^5 

X 

0 

1 

ro 

ro 

ro 

■-t=- 

.11 x 1o 2 .12 2.8 

.21* x io" 2 .27 1.6 

q 1 

*2 

.57 >< io -1 1.1 Y 

.57 x io " 1 1.1 f \ 

q l 

*6 

.12 x io" 2 .24 1.3 

.28 x io" 2 .35 1.7 

q 2 

q 6 

. 16 x io" 2 .2o 1 .9 

.28 x io" 2 .33 .98 

*2 

q 6 

.57 x Io" 1 It. It \ / 

.57 x Io" 1 4.1 Y 

.26 x io" U .55 /\ 


Finally we give an example in which the penalty method was used 
effectively. Again we consider (5- 1) and the "observations" are 
calculated from (5.2). In this case we take as initial guess 
(.3,0,0). This leads to a highly oscillatory behaviour of the 
trajectories of the differential equation and requires to 
decrease the tolerance in the Gear algorithm to io"'’. With 
this small value for the tolerance the three-level Levenberg- 
Marquardt became ineffective and in Table 6 we therefore show 
the numerical data that were obtained after the first level. 

The penalty function ip was taken to be ip = ip^ + ip 2 + > 


where 
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for 


q x * b 


M^i) = 


(q 1 -a) for q a < a 

(qj-b) 2 for q a > b 


with \p 2 and <j> 6 defined analoguously . The "correct value" 
for the parameter vector is = > as before 


Table 6. Optimal parameter and relative error i.p. for simultaneous 
(q^ ,q^ ,qg )-search with use of penalty functional (N=*t). 
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The numerical results of Table 6 were obtained with the scaling 
factor (modifying the Marquardt parameter) set initially at the 
same value (2.0) as in all other examples, in this case the 
scheme without penalization term does not converge. Convergence 
of the scheme without penalization can be obtained by taking 
the scaling factor as small as l.lj in this case the optimal 
parameters do not differ essentially from those of Table 6 when 
penalization is used. This and other examples demonstrate that 
by using a penalty functional in the fit to data criterion some 
robustness in the choice of the parameters in the Levenberg- 
Marquardt algorithm is gained. 
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